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Foreword 


The  Defense  Modeling  and  Simulations  Organization  (DMSO),  through  the  Air  Force 
Combat  Climatology  System  Modeling  and  Simulation  Division  (AFCCC/MS)  is  supporting  the 
Air  Force  Research  Laboratory  (Phillips  Lab)  contract  F19628-97-C-0025.  This  contract  is  a 
three-year  project  to  develop  a  method  for  estimating  climate  statistics  at  high  spatial  resolution 
in  both  data  sparse  and  data  rich  regions.  The  method,  given  the  name  CLIMOD,  is  based  on  the 
use  of  a  numerical-dynamical  atmospheric  model  to  generate  a  large  sample  of  simulated 
atmospheric  data  from  which  climatological  statistics  are  computed.  The  objectives  of  this 
project  are  to:  (1)  demonstrate  the  robustness  and  general  applicability  of  the  method  by 
generating  and  evaluating  climatological  databases  for  all  seasons  in  two  climatologically 
distinct  regions  of  the  world  that  are  of  interest  to  DOD  operations;  (2)  refine  the  CLIMOD 
technique  to  optimally  blend  available  observations  and  numerically  simulated  data  to  produce 
the  highest  quality  climate  statistics;  (3)  Extend  the  technique  to  provide  estimates  of  derived 
variables  and  additional  statistical  parameters  that  are  frequently  required  by  DOD  applications; 
and  (4)  construct  an  integrated  CLIMOD  software  package  that  can  be  efficiently  used  to 
generate  high  quality  simulated  climate  statistics  for  any  location  in  the  world. 

The  work  during  the  first  year  of  the  project  includes:  (1)  development  of  a  method  to 
optimally  blend  the  available  observations  and  model-generated  data  into  a  regional  statistical 
climatology;  (2)  development  of  a  method  to  objectively  assess  the  quality  of  the  simulated 
climate  statistics  in  data  rich  and  data  sparse  regions;  and  (3)  the  generation  of  three-dimensional 
simulated  climate  statistics  for  all  months  of  the  year  over  a  10-year  period  for  two  different 
regions  of  the  world:  Korea  and  the  Mideast. 

The  second  year  of  the  project  will  focus  on:  (1)  the  execution  of  data  denial  experiments 
to  assess  the  impact  of  each  data  type  on  the  quality  of  the  climate  statistics  generated  by  the 
CLIMOD  approach;  (2)  the  development  and  implementation  of  methods  to  compute  the 
statistics  of  derived  variables  from  the  basic  meteorological  variables  produced  by  the  model;  (3) 
the  use  of  the  derived-variable  computational  procedures  to  generate  climatological  datasets  for 
the  derived  variable  for  each  of  the  two  regions  for  which  a  climatology  of  basic  variables  was 
created  in  the  first  year;  and  (4)  the  development  of  methods  to  compute  specialized  statistics 
such  as  the  probability  of  extreme  events  from  the  basic  CLIMOD  output. 

The  work  during  the  third  year  will  consist  of:  (1)  the  documentation  and  delivery  of  the 
derived  variable  climatological  data  to  AFCCC;  (2)  the  construction  of  an  integrated  CLIMOD 
software  package  which  could  be  used  to  generate  climatological  databases  or  as  a  tool  to  further 
refine  the  CLIMOD  technique;  (3)  delivery  of  the  software  package  to  AFCCC;  and  (4)  the 
generation  of  three-dimensional  simulated  climate  statistics  for  all  months  of  the  year  over  a  10- 
year  period  for  a  tropical  region. 

This  interim  report  covers  Tasks  1  through  Tasks  4  of  this  project. 


1.0  Introduction 


There  are  multitudes  of  DoD  simulation  applications  that  require  statistical  values  of 
various  environmental  parameters  that  describe  the  local  3-D  climate  in  great  detail.  These 
applications  include  battle  to  global-scale  simulations,  real-world  contingency  ("what  if) 
situations,  ,  and  weapons  system  testing  and  development.  The  most  obvious  and  direct  way  to 
obtain  local  climate  statistics  is  to  calculate  them  from  long-term  point  observations.  However, 
the  usepfJong-term  observational  datasets  . imposes  at  least  two  limitations:  (1)  long-term  point 
observational  datasets  are  not  currently  available  for  many  areas  of  the  world;  and  (2)  the 
representativeness  and  quality  of  observations  change  in  time  as  observing  sites  are  relocated,  the 
land  use  around  a  site  is  changed  or  new  instrumentation  is  used.  An  alternative  is  to  generate 
statistic^  from  a  set  of  simulations  using  a  Mmited-area  high  resolution  mesoscale  model  that  is 
based  upon  the  physical  processes  of  the  atmosphere.  This  approach  can  be  used  to  address  two 
different  t3q)es  of  applications.  One  type  of  application  is  to  use  the  high  resolution  model  to 
infer  a  set  of  climate  statistics  from  a  climatic  state  that  has  been  generated  by  a  General 
Circulation  Model  (GCM).  The  high  resolution  model  is  used  to  add  additional  detail  to  the 
GCM  simulation.  In  this  approach  the  high  resolution  model  is  essentially  a  nested  model  within 
the  GCM.  The  initial  and  boundary  condition  data  are  extracted  from  the  GCM  output.  Giorgi 
et  al.  (1993)  used  a  version  of  the  Penn  State/NCAR  mesoscale  model  to  execute  a  2-year 
simulation  that  used  initial  and  lateral  boundary  data  extracted  from  a  GCM.  This  work  was 
recently  extended  by  Giorgi  and  Marinucci  (1996). 

In  the  second  type  of  application,  which  is  used  in  this  project,  the  objective  is  to  estimate 
the  actual  climate  statistics  for  a  particular  period  of  time  over  a  specified  region.  This  is  the 
objective  that  is  most  relevant  to  DQD  simulation  activities,  because  it  can  be  used  to  simulate 
actual  periods  of  records.  This  objective  can  be  addressed  by  executing  the  model  for  a  long 
period  pf  time  in  a  data  assimilation  mode  and  allowing  the  model  to  fiise.the  available  data  into 
a  3-D  climatological  dataset.  In  this  mode  the  model  continuously  or  periodically  ingests  the 
availably  observed  data  through  one  of  several  possible  data  assimilation  techniques  (Newtonian 
relaxation,  periodic  reanalysis,  etc.).  The  model  then  dynamically  fuses  the  available 
observations  from  the  region  of  interest  with  its  knowledge  of  the  surfece  characteristics  of  the 
earth  and  the  basic  principles  of  physics  to  generate  estimates  of  local  climate  statistics  at  all 
points  pn  a  three  dimensional  grid  that  covers  the  xegion.  This  technique  has  been  given  the 
name  CLImate  statistics  by  a  dynamical  MODel  (CLIMOD).  The  climate  results  are  likely  to  be 
much  more  representative  of  actual  periods  of  record  than  results  from  the  GCM-nested  method. 
The  CLIMOD  technique  is,  in  general,  a  more  difficult  approach  than  the  GCM-nested 
application  because  the  observed  data  must  be  inserted  into  the  simulation  in  a  way  in  which  it 
will  beneficially  impact  the  resulting  statistics.  This  is  not  a  trivial  problem  because  the  data  has 
quality  apd  representativeness  attributes  that  can  sometimes  cause  it  to.  conflict  with  the  values 
generated  by  the  model. 

The  Air  Force  Research  Laboratoiy  funded  project  ,  through  the  Air  Force  Office  of 
Scientific  Research  (AFOSR),  (AFOSR  contract  F4960-95-1-0523)  has  established  a  proof-of- 
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concept  for  the  CLIMOD  method.  This  project  performed  a  set  of  experiments,  which 
determined  that  the  climate  statistics  produced  by  the  CLIMOD  technique  are  of  very  good 
quality  if  properly  computed.  However,  the  research  also  showed  that  the  data  assimilation 
strategies,  PEL  schemes,  micro  physics  and  the  cumulus  parameterization  schemes  employed  in 
the  model  greatly  influenced  the  quality  of  the  statistics.  This  current  research  has  demonstrated 
significant  success  in  implementing  new  data  assimilation  strategies,  PEL  schemes  and  in 
selecting  the  appropriate  micro  physics  and  cumulus  parameterizations  for  producing  good 
quality  climate  statistics. 

The  numerical-dynamical  model  approach  researched  in  this  project  has  proven  to  have 
many  advantages  over  the  traditional  interpolation  approaches  of  generating  climatological 
statistics  from  available  observations  including: 

(1)  the  numerical-dynamical  method  has.the  capability  of  providing  reasonable  estimates  of 

climate  statistics  in  areas  without  any  direct  point  observations; 

(2)  global-scale  gridded  analysis  data  exists  for  the  past  forty  years  to  support  applications 
of  the  method; 

(3)  a  dynamically  consistent  3-D  climatology  can  be  produced  for  any  location  on  the 
earth; 

(4)  output  can  be  used  to  calculate  parameters  for  hydrology,  aerosols  and  E-0  weapon 
systems; 

(5)  it  is  a  flexible  tool  which  can  be  configured  to  provide  information  at  virtually  any 
required  spatial  and  temporal  scale; 

(6)  the  method  can  be  easily  integrated  with  other  databases  to  support  the  creation  of 
"what  if  scenarios;  and 

(7)  it  can  be  used  to  assess  the  impact  of  future  changes  in  the  environment  such  as 
urbanization,  changes  in  land  use  and  battlefield  conditions  on  the  climatological 
statistics. 


The  following  personnel  are  involved  in  the  research. 
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Dr  JohnZack 
Dr  Glenn  Van  Knowe 
Dr  Charles  Graves 
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Capt.  Matt  Doggett 
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AFRL  Contract  Monitor,  AFRLA^SEE 
Technical  Monitor,  AFCCC/SYT 
Technical  Monitor,  AFCCC/CCX 


2 


2.0.  DESCRIPTION  OF  THE  CLIMOD  APPROACH 

The  objective  of  the  CLIMOD  technique  is  to  estimate  the  actual  climate  statistics  for  a 
particular  period  of  time  over  3.  specified  region.  The  approach  used  executes  a  numerical  model 
for  a  long  period  of  time  in  a  data  assimilation  mode  and  allows  the  model  to  fuse  the  available 
data  into  a  3-D  climatological  dataset.  A  study  done  in  conjunction  with.this  research  has  shown 
that  a  recently  developed  form  of  data  assimilation  called  Incremental  Update  Analysis  (lAU) 
(Stauffer  et  al.,  1991)  produces  the  best  results  of  the  several  data  assimilation  schemes 
evaluated.  Through  the  lAU  process  combined  with  the  model's  knowledge  of  the  surface 
characteristics  of  the  earth  and  the  basic  principles  of  physics,  the  available  observations  of  any 
specified  region  are  dynamically  fused  to  generate  estimates  of  local  (point)  climate  statistics  at 
other  locations  in  their  region  for  which  no  observational  climate  data  is  available.  Once  the 
basic  climate  database  is  created,  observed  climatologies  are  used  to  improve  the  mean  values 
and  variability  by  compositing  the  observed  data  with  the  simulated  data. 

2.1  GENERAL  SYSTEM  DESIGN. 

The  CLIMOD  method  is  designed  to  produce  the  simulated  climatologies  in  a  six  step 
process,  (1)  The  MASS  model  is  used  to.  simulate  10  years  of  hourly  observations  at  a  10  and 
40  km  grid  resolution.  The  simulations  are  generated  in  data  assimilation  mode  with  a  nested 
grid  system.  A  coarse  grid  with  a  3-D  grid  point  matrix  size  of  80x80  horizontal  grid  points  on 
20  vertical  levels  uses  a  40  km  horizontal  grid  spacing,  while  the  fine  mesh  grid  uses  the  same 
3-D  grid  point  matrix  size  but  employs  a  horizontal. grid  increment  of  10  km.  Thus,  the  coarse 
mesh  grid  covers  an  area  of  approximately  3200  by  3200  km  while  the  fine  mesh  operates  on  a 
theater  scale  area  of  approximately  800  by  800  km.  Attachment  1  outlines  the  steps  performed 
in  the  running  of  the  model.  (2)  The  model  output  is  qualitatively  checked  through  the  runs  for 
meteorological  soundness.  (3)  Observational  data  is  then  blended  with  the  simulated  data  to 
create  a  seamless  database  that  maintains  the  best  quality  at  every  location.  (4)  Climate  statistics 
are  then  generated  for  grids  and  specific  locations  within  the  fine  mesh  and  coarse  mesh 
domains.  (5)  Rigorous  quantitative  evaluation  is  made  by  comparing  the  composited  climate 
statistics  with  observed  climate  data.  This  is  done  using  a  method  developed  for  this. project  to 
establish  the  confidence  the  user  can  have  in  the  data.  (6)  The  data  is  then  put  into  a  standard 
format  for^storage. 

The  global  grid  point  analysis  dataset  (e.g.  the  NCAR/NCEP  GDAS  Reanalysis  dataset) 
and  all  available  rawinsonde  and  surface  data  are  being  used  to  generate  the  climate  grids.  It  is 
also  theoretically  possible  to  assimilate  other  types  of  data  such  as  radiance  data  or  derived 
sounding  m  cloud  information  fi:om  satellite  sensors.  The  use  of  such  data  may  have  a 
significant  beneficial  impact  on  certain  components  (e.g.  cloud  climatology)  of  the 
climatological  datasets.  However,  the  use  of  such  datasets  is  beyond  the  objectives  and 
resources  of  this  project.  Nevertheless,  this  does  represent  an  area  in  which  future  enhancements 
can  be  pade  to  the  CLIMOD  method. 
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2.2  DATA  DESCRIPTION. 


a.  Databases  Defining  the  Surface  Characteristics 

Land  Water  Distribution.  One  of  the  most  basic  grid  parameters  that  needs  to  be 
determined  is  whether  a  grid  box  is  composed  of  land  or  water.  An  integer  indicator  or  flag  is 
used  to /ieteamine  land  or  water  where  "0"4ndiGates  water  and  "1"  indicates  land.  The  datab"ase 
used  to  determine  the  land/water  flag  is  a  U.S.  Navy  dataset  obtained  from  the  National  Center 
for  Atmospheric  Research  (NCAR).  The  resolution  of  the  dataset  is  approximately  18  km  in  the 
mid-latitudes. 

Terrain  Elevation.  The  database  used  to  describe  the  terrain  elevation  for  each  model 
grid  box  is  the  U.S.  Central  Intelligence  Agency  5  minute. global  terrain  dataset.  The  resolution 
of  this  data  base  is  approximately  9  km. 

Sea  Surface  Temperature  (SST).  The  MASS  uses  SST  data  from  global  sea  surface 
temperature  climatological  datasets.  One  dataset  is  used  for  each  month  of  the  year  at  a  1  degree 
latitude-longitude  resolution. 

Land  use  and  Land  cover.  The  MASS  land  use/land  cover  dataset  is  based  on  the 
Anderson  Level  II  classification  scheme  (Anderson  et.  al.,  1976)  developed  by  the  U.S. 
Geological  Survey  (USGS).  The  Anderson  Level  II  land  use  scheme  provides  10  basic  land  use 
choices,  such  as  urban,  agricultural,  range,  wetland,  etc.  and  further  divides  each  major  category 
to  give  a  total  of  approximately  50  choices.  The  resolution  of  the  land  use/land  cover  is  5  x  5 
degree  latjtude-longitude  sectors. 

Soil-t>T)e  Classification.  A  soil  classification  database  is  also  used  in  the  MASS. 

formalized  Difference  Vegetation  Index  (NVDI).  Data  providing  information  on  the 
vegetation  is  also  used  in  the  MASS.  The  basic  data  has  been  collected  by  the  Advanced  Very 
High-Resolution  Radiometer  (AVHRR)  instrument  The  NVDI  used  is  derived  by  calculating  the 
normalized  difference  between  two  AVHRR  channels,  the  near  infrared  and  the  visible.^he 
resolution  of  ND  VI  dataset  is  5x5  degree  latitude-longitude  sectors. 

b.  Atmospheric  Databases 

Research  was  done  on  the  potential  of  using  the  NCEP/NCAR  40  Year  Reanalysis 
Project  datasets  _ from  NCAR  as  an  improved  source  of  data  to  initialize  and  provide  lateral 
boundary  conditions.  The  Reanalysis  Project  data  is  available  every  six  hours  versus  every  12 
hours  for  Global  Optimum  Interpolation  (GOI)  data  datasets  used  in  the  proof-of-concept  grant. 
MESO  had  to  first  develop  a  reader  that  converts  the  Reanalysis  Project  data  format  into  a 
MASS  readable  format.  MESO  then  performed  tests  on  the  Reanalysis  Project  data  comparing 
its  quality  with  the  previously  used  Global  Optimum  Interpolation  datasets.  It  was  found  that 
using  the  Reanalysis  data  better, limits  climate  drift  of  the  model  and  is  an  improvement  over  the 
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GOI  datasets.  The  specific  dataset  used  in  the  project  is;  ds090.0  1  January  1973  to  31  Dec 
1982  from  NCEP/NCAR  Reanalysis  Archive. 

The  10  years  of  rawinsonde  data  were  also  obtained  from  the  NCAR  dataset:  ds353.4 
RAOBS  ^PIBALS  Global.  Data  for  01  Jan  1973  -  31  Dec  1982. 

The  surface  data  used  was  obtained  from  the  Air  Force  Combat  Climatology  Center 
(AFCCC)  DATSAV  data.  A  reader  was  constructed  to  convert  the  DATSAV  into  a  format 
readable  by  MASS. 


2.3  DESCRIPTION  OF  THE  MESOSCALE  MODEL 

The  simulation  model  used  in  this  project  is  version  5. 1 1 .2  of  the  Mesoscale  Atmospheric 
Simulation  System  (MASS)  model  (Kaplan  et  al.  1982;  MESO,  1995;  Manobianco  et  al.,  1996). 
The  MASS  is  a  three-dimensional  primitive  equation  hydrostatic  numerical  model. 

a.  Basic  Model  Equations  and  Grid  System 

The  basic  model  formulation,  equations  and  grid  system  are  fully  described  in  the  MASS 
Version  5.6  Reference  Manual  (MESO,  Inc.,  1994).  The  MASS  uses  the  vertical  coordinate 
described  as 

ai) 

p.-p, 

where  p  is  the  pressure  at  the  level,  ps  represents  pressure  at  the  surface  and  pt  represents  pressure  at 
top  of  model  domain.  The  MASS  can  be  run  using  either  a  polar  stereographic  or  mercator  map 
projection.  For  this  experiment  the  polar  stereographic  projection  was  used  employing  the  map 
scale  faptor  m  =1  +sin((l))/I  +  sin((t>o),  where  (j)  .  denotes  any  latitude  and  <|)o  denotes  the  reference 
latitude.  The  MASS  uses  an  unstaggered  horizontal  grid  and  a  staggered  vertical  grid.  The 
atmospheric  portion  of  the  MASS  consists  of  a  hydrostatic  system  of  seven  model  prognostic 
variables;  u-wind  component,  v-  wind  component,  pressure,  temperature,  humidity  as  water  vapor 
mixing  ratio,  precipitation  moisture,  and  cloud  moisture. 

Horizontal  finite  differencing  is  done  using  a  fourth  order  scheme.  Time  marching  is 
accomplished  by  using  the  split  explicit  scheme.  A  forward  backwards  scheme  is  used  for  the 
capture  of  inertial  gravity  waves.  Advection  is  handled  using  the  Adams-Bashford  scheme 
(Meisinger  and  Arakawa,  1976). 


The  MASS  surface  energy  budget  is  formulated  after  Noilhan  and  Planton  (1989),  which  is 
based  on  the  original  work  of  Blackadar  (1976,  1979).  The  MASS  calculates  the  major  impacts  of 
both  short  and  longwave  radiative  transfer  on  heating  the  ground  surface,  clouds  and  atmosphere. 
There  afe  three  radiative  situations  (clear,  cloudy  ,and  partly  cloudy)  considered  in  the  MASS  for 
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long  and  shortwave  radiation.  The  fractional  sky  cover  is  estimated  based  on  a  relative  humidity- 
cloud  fraction  relationship  at  each  model  level.  Walcek  (1994)  calculated  the  observed  cloud  cover- 
relative  humidity  relationships  for  an  early  spring  midlatitude  cyclone  over  land.  The  main 
observed  difference  is  the  approximate  20%  increased  cloud  cover  with  relative  humidities  between 


The  solar  radiation  wavelengths,  referred  to  as  shortwave  radiation,  considered  in  the  MASS 
ranges  from  0.3-4  micrometers.  The  solar  radiation  is  further  divided  into  a  visible  (0  3-0  75 
micrometer)  and  an  infrared  (0.75-4.0)  micrometer  band.  The  teirestrial  or  longwave  radiation 

ranges  from  4.0  to  50  micrometers  and  accounts  for  longwave  radiative  transfer  from  both  the 
ground  and  the  clouds. 

For  shortwave  absorption  in  clear  skies  the  direct  radiative  heating  of  the  atmosphere 
(shoi^ave  radiative  flux  divergence)  is  ignored.  Only  the  radiative  heating  due  to  the  shortwave 
irradiance  received  at  the  ground  is  parameterized.  Rayleigh  scattering  and  Mie  scattering  is 
accounted  for  at  the  rnolecular  and  particulate  size  level  as  the  shortwave  radiation  passes  through 
the  atmosphere.  The  impact  on  the  zenith  angle  due  to  terrain  slope  is  also  taken  into  account.  The 
basic  relationship,  integrated  at  each  model  level,  used  to  calculate  shortwave  radiation  reaching  the 
ground  during  clear  sky  conditions  is  Rsfrk)  =  (1  +  At  -Aa)Rsfrk+l);  where  R^frk)  denotes  radiation 
reaching  level  k  from  the  level  above  and  k+1  At  represents  the  change  in  transmissivity  from  k  to 
defines  the  change  in  absorptivity  from  k  to  k+1.  For  cloudy  sky  conditions,  the  shortwave 
radiatiop  is  partitioned,  into  the  visible  and  infrared  bands.  The  total  solar  radiation  reaching  the 
ground  during  cloudy  conditions  is  calculated  by  Rsfrcloudy)  =  tvisRsofrvis)  +  ti,R,ofrir),  where  fris 
IS  the  tr^missivity  m  the  visual  range  of  the  radiation  spectrum,  R«,frvis)  represent  the  solar 
constant  at  the  top  of  the  atmosphere  in  the  visual  spectrum,  tir  represent  the  transmiwity  in  the 
infrared  spectrum,  and  Rsol(ir)  represent  the  solar  constant  at  the  top  of  the  .  atmosphere  in  the 
infrared  part  of  the  spectrum,.  For  partly  cloudy  conditions,  the  clear  and  cloudy  downward 
shortwave  irradiances  are  weighted  by  the  total  column  cloud  fraction  as  follows: 

RsJ,  (partly  cloudy)  =  (1-fraCooi)  Rsoi  (clear)  +  fraCcoiRsoi(cloudy), 
where  Rso  is  the  solar  constant,  fraCcoi  equals  cloud  fraction  for  the  entire  column 

The  parameterization  of  longwave  radiative  transfer,  used  to  account  for  longwave  cooling 
and  heating  of  the  ground  and  atmosphere,  uses  the  broadband  approach  of  Sasamori  (1972)  as 
described  by  Pielke  (1984).  The  effects  of  water  vapor  and  carbon  dioxide  longwave  radiation 
absorptipn  are  parameterized.  In  the  cloudy  sky  conditions  the  effective  emissivity  of  the  cloud  is 
also  taken  into  account.  In  the  partly  cloudy  situation,  the  clear  and  cloudy  longwave  fluxes  are 
combined.  The  longwave  radiative  process  is  partitioned  into  downward  ARa  and  upward 
ARLtfluxes  for  the  clear  and  cloudy  cases. 


b.  MASS  Planetary  Boundary  Layer 

Two  options  exist  for  the  MASS  model  planetary  boundary  layer  formulation.  The  first 
closely  follows  the  high  resolution  PBL  model  described  by  Zhang  and  Anthes  (1982),  which  is 


mostly  based  on  the  work  of  Blackadar  (1976,  1979).  This  scheme  has  been  used  for  many  years 
by  several  modeling  groups  with  good  results.  Some  of  the  changes  in  the  scheme  from  the 
original  Zhang  and  Anthes  paper  are  derived  from  changes  made,  in  the  MM4  version  of  the  Penn 
State/NCAR  mesoscale  model  (Stauffer,  1992).  A  vertical  profile  of  temperature,  water  vapor 
mixing  ratio,  and  u-  and  v-components  of  wind  on  an  arbitrary  set  of  pressure  levels  are  passed 
into  the  scheme,  along  with  relevant  surface  characteristics.  The  result  of  the  scheme  is  the 
calculation  of  time  tendencies  of  each  of  the  four  variables  listed  above,  along  with  some 
diagnostic  PBL  parameters,  such  as  the  depth  of  the  boundary  layer. 

The  second  option  is  based  on  the  work  of  Therry  and  LaCarrere  (1983).  It  is  a  so-called 
1.5  order  closure  scheme  in  which  there  is  a  prognostic  equation  for  turbulent  kinetic  energy 
(TKE)  and  the  turbulent  mixing  coeffiecient  is  a  function  of  TKE.  The  boundary  conditions 
were  taken  from  Mailhot  and  Benoit  (1982). 

The  MASS  model  divides  the  Planetary  Boundary  Layer  (PBL)  into  two  layers.  The 
lower  layer  represents  the  surface  or  constant  flux  layer,  the  upper  or  above  surface  layer 
represents  the  Ekman  layer.  The  Blackadar  and  TKE  schemes  are  virtually  identical  in  the 
surface  layer.  The  differences  are  in  the  above  surface  layer.  The  MASS  surface  layer 
formulation  is  founded  on  separating  the  PBL  into  four  regimes  based  on  stability  to  calculate 
the  scaling  parameters.  The  criteria  used  to  choose  among  three  stable  regimes  and  one  unstable 
regime  is  the  ratio  of  the  PBL  height  (h)  and  the  Monin-Obukhov  length  (L)  combined  with  the 
bulk  Rich^dson  number  (Rb). 

The  critical  value  of  h/L  =  1.5  is  used  to  separate  unstable  (free  convective)  from  three 
stable  regimes.  If  h/L  is  equal  to  or  less  than  1.5  than  Rb  is  used  to  identify  which  stable  regime  is 
used.  The  totally  stable  regime  occurs  when  Rb  is  greater  than  0.2.  The  damped  mechanical 
turbulence  regime  is  used  when  Rb  is  between  0  and  0J2.  When  Rb  is  equal  to  or  less  than  zero,  a 
forced  convective  regime  is  invoked.  Each  regime  has  different  values  of  heat,  moisture  and 
momentum  coefficients  as  described  in  the  MASS  Reference  Manual  (1995). 

The  MASS  model's  three  stable  and  one  unstable  (convective)  surface  PBL  regimes  are 
calculated  using  the  following  method.  First,  profiles  of  four  variables  are  used  as  input:  tempera¬ 
ture,  moisture  (mixing  ratio),  u  wind  component  and  the  v  wind  component.  The  criteria  for 


Rb=(g  Za^^a)(ea-eg)/(Vaf  <  (2.2) 

selection  of  n.  regime  is  based  on  Rb  and  h/L,  Specifically, 
where:  Za  equals  height  of  first  model  level 

g  is  acceleration  of  gravity 

0a  equals  potential  temperature  at  the  first  model  level 
0g  equals  potential  temperature  at  the  equals  ground  surface  and 
Va  equals  wind  velocity  at  the  first  model  level. 

For  the  three  stable  regimes  the  scaling  parameters  for  heat,  moisture  and  momentum  are  set 
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to  be  equal.  That  is  x|/h  -  Vq  =  \|/m,  where  h  denotes  heat,  q  denotes  moisture,  and  m  denotes 
momen^m.  In  the  MASS  «nly  and  vj/m  are  used.  The  Rb  criteria  distinguishes  which  stable 
layer  is  used^  Two  stable  basic  types  of  stable  layers  are  possible.  If  Rb  >  0,  stability  is  indicated 
with  little  wind  shear.  If  Rb  >  0.^2,  the  most  stable  regime  results,  vinth  only  weak  turbulence 
allowed.  For  the  most  stable  regime,  the  frictional  velocity  u*  =  0.1  (land)  and  0.0  (water)  and  the 
scaling  parameter  hi  =  -10  ln(z3/zo).  If  0  >  Rb  <  0.2  (between  0  and  0.2),  then  the  damped 
mechanic^  turbulent  regime  results  where  the  scaling  parameter  x]/  =  -5(Rb/(1.10-5Rb)ln(za/2o). 

If  Rb  <  0,  you  must  go  to  h/L  to  determine  if  it  is  a  stable  (forced  convective)  or  unstable 
(free  copvective)  regime.  If.h/L  is  less  than  1.5,  forced  convection  results.  With  forced,  convection, 
the  scaling  parameters  are  set  equal  to  0  and  the  boundary  layer  depth  is  set  to  the  first  level  above 

the  surface  layer  where  0  >  0a.  Ifh/L  is  greater  than  1.5,  unstable  free  convection  results  with  the 
scaling  parameters  not  set  to  be  equal. 


The  stable  regime  tendencies  are  calculated  from  the  following  basic  relationship: 


dx  _(X^-KmX,) 
dt  pA  zi 


(2.3) 


where  x,  represents  any  variable,  Xgis  the  value  of  the  variable  at  the  surface,  xi  is  the  value  at  the 
first  model  level,  K*  is  the  eddy  viscosity  defined  in  "K-Theory"  (Stull,  1988X  p  is  the  density,  and 

Azl  is  the  depth  of  the  first  model  level.  Heat  flux  at  the  top  of  the  layer  is  calculated  from  a  linear 
approach. 


For  the  unstable  free  convective  regime  the  tendency  equations  are  based  on  the  principle 
that  a  fraction  of  the  total  mass  is  being  transferred  between  the  surface  layer  and  each  PEL  layer. 
The  fraction  of  mass  (m)  exchanged  between  any  level  and  the  surface  layer  per  unit  time  is 
as  follows: 


'  Zm 

PCp(l-  Em)  \(Sa-  6(z))dz 


where  is  an  entrainment  parameter  set  to  0.2  and  Zm  is  the  height  of  the  level  where  a  rising 
parcel  of  air  loses  its  buoyancy.  The  heat  flux  is  calculated  using  a  non-linear  term. 

The  Blackadar  and  TKE  scheme  have  significantly  different  formulations  in  the  above 
surface  or  Ekman  layer.  The  Blackadar  sceme  uses  the  K-theory  approach  as  described  by 
Black^^  (1976).  All  of  the  stable  layers  use  the  K-theory  approach  in  the  Ekman  layer  where  Km 
—  Kh  —  Kq  is  assumed  which  relates  the  turbulent  flux  to  the  associated  mean  variables.  In  the 
MASS  the  vertical  diffusivity  coefficient  terms  Kh  and  Km  are  calculated,  from  use  of  the  second 
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order  she^  term  and  the  local  Richardson  number.  Specifically,  Km  is  calculated  from 


Xzo  +  Ri  <  Ric 

Km  =  {  Ric 

K:0,if  Ri^  Ric  ■  (2.5) 

Where  Ri  is  the  local  Richardson  number,  Ric  is  the  critical  Richardson  number,  Kzo  is  the 
background  diffiisivity,  S  is  the  vertical  wind  shear,  1  is  the  mixing  length  set  to  100  m.,  k  is  the  von 
Karman  constant.  The  value  of  X  is  obtained  from  the  K-Theory  relationship  which  is  u'X  =  -K 

dX/dz  where  X  denotes  any  of  the  variables. 

For  the  free  convective  regime  the  local  exchange  coefficient  is  used  above  the  mixed  layer. 
Within  the  mixed  layer  the  tendency  equations  are  based  on  the  fractional  mass  (m)  exchanged 
between  layers  described  in 


dX/dt  =  m(Xa-X}, 


(2.6) 


where  X  rejiresents  any  of  the  model  variables. 


TKE  scheme  calculates  the  tendencies  above  the  surface  layer  using  the  prognostic  equations 
for  the  l.Sorder  closure  scheme  as  outlined  by  Stull,  (1988). 


(2.7) 


dw 


5(v'W) 

Sz 


ae  ,  ff  ae  =  i  j  _  LvE  _  a^v) 

~St  pCp  axj  pCp  az 


ac  ,  rr  ac  ^  Sc  _  a^’wQ 

"ar  ^  j  Pair  az 


(2.8) 


<2.9) 


(2.10) 
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(2.11) 


de 


+  U 


^Uj  +e)] 

~Bz  dz 


The  general  notation  is  as  follows: 


1. 


2. 


4. 


Mean  quantities  are  generally  capitalized  with  an  overbar  (-).  They  are  the 
west-to-j^ast,  southTto^north  and  vertical  components  of  velocity  (U,  V, 
potential  temperature  0,  pressure  P,  tracer/moisture  variable  mixing  ratio  c, 
density (p  and  turbulent  kinetic  energy  e. 


Perturbations  from  the  mean  for  a  given  volume  for  the  above  variables  are 
denoted  by  a  prime  ( ' ),  for  example,  u'; 

The  yalue  of  the  perturbation  quantities  derived  through  a  Reynolds  averaging 
procedure  have  an  overbar,  for  example  u'  or  u'w; 

Some  terms  use  a  subscript  j  to  represent  all  three  cartesian  directions.  For 
example,  Uj  refers  to  the  three  wind  components  and  xj  refers  to  the  three 

directions.  Consequently,  Uj^.is  shorthand  for 


The  terms  on  the  left  hand  side  in  each  equation  represent  the  local  tendency  and  advection, 
respectively.  Equations  (2.7)  and  (2.8)  describe  the  evolution  of  the  two  horizontal  components 
of  momentum.  The  first  term  oathe  right  hand  side  represents  the  effects  of  the  Coriolis  and 
pressure  gradient  forces.  The  last  term  represents  turbulent  fluxes.  Equation  (2.9)  describes  the 
evolution  of  the  potential  temperature.  The  terms  on  the  right  hand  side  represent  the 
contributions  of  radiation,  water  phase  changes  and  turbulent  fluxes,  respectively.  Equation 
(2.10)  dpcribes  the  evolution  of  a  particular  phase  of  water  (vapor,  solid,  liquid)  or  some  other 
passive  tracer  C.  The  first  term  on  the  right  hand  side  represents  sources  or  sinks  (condensation, 
evaporapon,  etc)  while  the  second  represents  turbulent  fluxes.  Equation  (2. 1 1)  is  the  TKE 
budget  equation.  The  terms  on  the  right  hand  side  represent  buoyant  production/loss,  mechanical 
or  shearjjroduction/loss,  turbulent  transport  and  redistribution  by  pressure  perturbations,  and 
viscous  dissipation,  respectively. 

The  general  expression  for  the  parameterization  of  the  second  order  fluxes  and  other 
terms  are; 


u'w'=-K 


8U 


(2.12) 


vw 


-K 


(2.13) 


c'w'.-Kh^ 


(2.14) 


w' 


(2.15) 


The  second  term  on  the  right  hand  side  of  (10-97)  allows  for  counter  gradient  buoyancy  fluxes 
under  unstable  conditions. 

The  eddy  diffusivity  coefficients  and  dissipation  of  TKE  are; 

Kn,  =  0.51^1/2  ,Kh  =  1.35K,  (2.16) 

a3/2 

and 

lie  and  Ig  Ik  and  Ig  Therry  and  Lacarrere  (1983)  who  used  the  results  of  a  third  order  closure 
simulation  to  parameterize  them  as  variables  rather  than  constants.  They  parameterize  Ik  and.  Ig 
as  as  rather  complicated  functions  of  the  PBL  depth,  the  local  height  above  ground  level,  the 
local  stability; 


ls  = 


(2.18) 


mi  = 


1 


1  +  0.005^ 
Kz 


m2  = 


1  _  ^mo 

h 


rimo’<  0 


(2.19a  and  b) 


OJ 


surface  layer  heat  flux  is  upward 
surface  layer  heat  flux  is  downward 


(2.20a  and  h) 


lk  = 


(2.21) 


1  0.0025-^ 

Jez 


mo  =  mo 


(2.22a  andb) 
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where  ?  is  the  height  above  the  ground  and  k  is  von  Karman's  constant  (0.4).  The 
countergradient  term  in  the  vertical  temperature  flux  term,  yc,  is  parameterized  as; 


Owl 


(2.23) 


w'9', 

where  0»  =  -W  and  w* 


W9 


^  'gw^’^1 

..  e  J 


1/3 


The  subscript  s  refers  to  the  surface  layer  while  z  refers 


to  height  above  ground  level.  Finally,  the  fourth  term  in  the  TKE  budget  equation,  is 
parameterized  as; 


iR-e  =  0.5|[l,e-g  -  V-wV] 


In  the  lowpt  model  layer  which  is  considered  to  be  the  surface  layer,  the  TKE  is  set  to; 

f  3.75ui  (stable) 

e  =  <  r-z\2r.  (2.25) 

®  ^.75ui  +  0.2wi  +  ui  (unstable) 

The  atmosphere  is  considered  to  be  unstable  if  potential  temperature  decreases  with  height 
between  (he  surface  and  the  lowest  model  layer  above  the  surface. 

c.  MASS  Soil  Hydrology 

The  soil  hydrology  scheme  is  based  on  Mahrt  and  Pan  (1984)  model.  The  soil  is  divided 
into  a  shallow  5  cm  deep  surface  layer  and  a  deeper  5-30  cm  layer.  The  plant  roots  are  assumed  to 
be  in  the  deeper  layer.  A  surface  "cover"  moisture  reservoir  is  parameterized  to  account  for 
intercept^  rainfall.  The  hydrology  scheme  uses  the  land  use,  soil  type  and  vegetation  index  data 
for  its  computations. 

The  soil  moisture  budget  takes  into  account  the  precipitation  that  falls,  the  amount  that  is 
intercep(ed.by  the  surface,  and  the  amount  that  infiltrates  into  the  soil.  The  model  also  accounts  for 
the  diffusion  of  moisture  in  the  soil,  and  the  evaporation  of  bare  soil  and  vegetation. 

The  effects  of  snow  cover  are  also  parameterized  in  the  model.  The  albedo  (which  is  by  far 
the  most  important),  evaporation  rates,  and  roughness  length  are  all  altered  to  account  for  snow 
cover. 

d.  MASS  Moisture  Physics 

The  MASS  handles  the  moisture  physics  and  latent  heating  on  two  scales.  The  larger  scale 
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is  handled  as  grid-scale  processes  explicitly  calculating  the  latent  heating  at  each  grid  point.  The 
convective  latent  heating  is  handled  through  sub-grid  scale  parameterization. 

The  grid  scale  moisture  physics  can  be  handled  in  two  ways.  One  method  employs  a 
diagnostic  scheme  precipitating  all  liquid  and  frozen  water  out  of  the  atmosphere.  This  is  the 
scheme  employed  in  the  40  km  horizontal  courser  mesh  nest.  The  other  ,  option  is  a  prognostic 
scheme,  which  accounts  for  the  phase  change  processes  involving  both  the  cloud  size,  water 
droplets,  and  ice  crystals  as  well  as  the  precipitation  size  raindrops  and  snow.  The  major 
disadvantage  of  using  the  prognostic  scheme  is  the  increased  computational  time.  The  prognostic 
scheme  is  used  in  the  10  km  horizontal  finer  mesh  nest. 

Convective  latent  heating  is  handled  through  the  use  the  of  Kuo-Anthes  (Kuo  1965;  Anthes 
1977)  parameterization  scheme.  The  Kuo-Anthes  scheme  uses,  the  principle  that  the  magnitude  of 
the  convective  scale  heating  is  proportional  to  the  grid  scale  moisture  supply. 

2.4  MASS  Configuration  for  CLIMOD. 

Several  issues  concerning  the  quality  of  the  data,  model  run  time  and  model  configuration 
were  determined  through  evaluating  the  results  of  the  proof-of-concept  project  as  they  pertain  to 
the  current  research.  The  primary  issues  investigated  and  resolved  are:  (1)  which  planetary 
boundary  layer  scheme  to  use,  (2)  which  cumulus  scheme  to  use,  (3)  moisture,  physics  issues  (4) 
and  how  to  solve  problems  during  periods  of  data  assimilation.  Much  of  the  early  testing 
focused  on  determining  the  “best”  model  configuration  for  the  chmate  runs.  Utilizing  the  test 
results,  modifications  were  made  to  the  model  making  it  more  optimal  for  climate  purposes. 

Because  the  focus  of  the  climate  project  is  to  improve  the  capability  of  simulating 
climatologies  near  the  earth's  surface,  a  comparison  was  made  using  the  two  MASS  PBL 
schemes;  the  traditional  Blackadar  scheme  and  the  TKE  PBL  scheme  based  on  turbulent  kinetic 
energy.  The  comparison  used  the  results  from  the  proof-of-concept  research  as  well  as  results 
from  additional  investigation  of  the  subject  as  part  of  this  project.  The  comparison  determined 
that  TKE  produced  superior  results  for  the  temperature,  dew  points  and  wind  speeds. 

It  was  observed  that  a  small  but  significant  drying  of  the  simulated  atmosphere  occurred 
each  time  rawinsonde  data  was  assimilated  into  a  model  run  through  the  Newtonian  Relaxation 
or  'nud^g"  scheme  based  upon  an  algorithm  that  is  similar  to  that  used  by  Stauffer  et  al.  (1991). 
This  resulted  in  a  strong  low-cloudiness  and  precipitation  bias  for  the  first  three  hours  of  the 
simulation  (0-3  UTC).  The  reason  for  this  was  determined  to  be  the  slight  low  humidity  bias  in 
the  rawinsonde  data  combined  with  the  nudging  method  over-forcing  the  model  solution  to  the 
dry  biased  observed  rawinsonde  data.  The  specific  details  of  the  simulation  strategy  that  is  being 
employed  to  generate  the  40  km  and  10  km  databases  were  determined  from  the  results  of  the 
analysis  of  simulation  experiments  that  were  completed  in  the  proof-of-concept.  The  issues 
investigated  were  (1)  the  types  of  data  that  will  be  assimilated  into  the  simulations;  (2)  the 
frequency  of  the  data  assimilation;  and  (3).  the  method  used  to  assimilate  the  data. 
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The  first  major  issue  of  data  assimilation  addressed  was  what  t5^e  of  data  to  use.  The 
global  grid  point  analysis  dataset  (e.g.  the  NCAR/NCEP  GDAS  Reanalysis  dataset)  was 
determined  to  be  the  best  available  choice  for  initialization  and  lateral  boundary  conditions.  The 
next  task  was  to  determine  which  was  the  best  strategy  of  using  the  available  rawinsonde  and 
surface  data.  There  were  some  issues  concerning  the  impact  of  data  assimilation  within  the 
planetary  boundary  layer.  For  example,  based  on  their  experiments  with  the  PSU/NCAR  model, 
Stauffet  .et  al.  (1991)  recommend  that  surface  temperature  data  not  be  used  in  a  Newtonian 
relaxation  data  assimilation  scheme  within  the  planetary  boundary  layer.  They  note  that  the 
continuous  assimilation  of  surface  temperature,  data  can  have  a  negative  impact  on  the  model's 
boundary  layer  parameterization  scheme.  Analysis  of  the  experiments  performed  in  the  proof- 
of-conc^pt  project  suggested  that  the  assimilation  of  surface  data  can,  at  times,  have  a  negative 
impact  on  the  near-surface  component  of  the  model  simulation.  However,  even  with  the 
problems  of  assimilating  observed  data,  it  was  found  that  the  addition  of  observed  surface  and 
rawinsonde  data  helps  remove  the  model  drift.  It  was  further  determined  that  it  is  theoretically 
possible  lo  assimilate  other  types  of  data  such  as  radiance  data,  derived  soundings,  or  cloud 
information  from  satellite  sensors.  The  use  of  such  data  may  have  a  significant  beneficial  impact 
on  certain  components  (e.g.  cloud  climatology)  of  the  climatological  datasets.  However,  the  use 
of  such  datasets  is  beyond  the  objectives  and  resources  of  this  project.  Nevertheless,  this  does 
represent  an  area  in  which  future  enhancements  can  be  made  to  the  CLIMOD  method. 

The  second  major  issue  addressed  was  the  frequency  of  data  assimilation.  An  intuitive 
approaph  would  be  to  assimilate  the  data  at  the  same  frequency  at  which  it  is  observed.  In  this 
approach,  rawinsonde  data  would  be  assimilated  at  12-hour  intervals  and  surface  data  would  be 
assimilated  at  hourly  intervals  in  regions  where  it  is  available  at  that  frequency,  and  at  3  hour  or 
6  hour  intervals,  in  other  areas.  However,  it  was  found  that  frequent  data  assimilation  can  have 
some  drawbacks  due  to  characteristic  problems  with  certain  types  of  data.  For  example, 
rawinsonde  data  tends  to  have  a  dry  bias  in  atmospheric  environments  that  are  near  saturation. 
Thus,  the  assimilation  of  rawinsonde  moisture  data  can  have  the  effect,  of  erroneously  reducing 
the  amount  of  cloud  cover  and  precipitation  produced  in  a  model  simulation.  Therefore,  the 
assimilafion  of  rawinsonde  moisture  data  at  more  frequent  intervals  can  generate  more  of  a  dry 
bias  in  a  simulation.  This  is  less  of  a  problem  in  forecast  simulations  because  the  model  will 
spin-up  to  the  proper  degree  of  saturation  in  a  few  hours  after  ,  initialization  and  then  be 
unencumbered  by  the  rawinsonde  dry  bias.  However,  in  a  simulation  to  generate  climatological 
datasets,  the  repeated  assimilation  of  rawinsonde  data  can,  for  example,  introduce  a  dry  bias  into 
the  cloud  and  precipitation  climatologies 

In  an  attempt  to  maximize  the  positive  aspects  of  using  observed  data  to  remove  model 
drift  and  to  minimize  the  negative  impact  of  observed  biases,  it  was  determined  that  the  impact 
of  different  data  assimilation  strategies  on  the  resulting  climatological  databases  needed  to  be 
investigated  further  as  part  of  this  project.  The  results  of  the  data  assimilation  experiments  were 
used  as  a  guide  for  selecting  the  strategy  to  be  used  for  the  generation  of  the  climatological 
databases  jfor  Korea  and  the  Middle  East. 

There  were  many  possible  candidates  for  the  four-dimensional  data  assimilation  method. 
However,  5ome  are  too  computationally  intensive  (e.g.  the  adjoint  technique)  or,  as  yet. 
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unproven  to  be  used  in  this  project.  The  two  methods  that  have  been  widely  used  and  are 
computationally  manageable  are;  (1.)  discontinuous  data  assimilation  through  periodic  objective 
analyses  using  the  simulation  output  as  a  first-guess;  and  (2)  continuous  data  assimilation 
through  the  use  of  Newtonian  relaxation  (also  known  as  "nudging").  The  discontinuous 
approach  has  been  widely  used  in  operational  forecast  centers  to  initialize  forecast  simulations. 
For  example,  the  12-hour  forecast  from  a  previous  run  is  used  as  a  first-guess  for  the 
initialization  analysis  for  the  next  forecast  simulation.  The  advantage  of  this  approach  is  that 
model  simulation  is  interrupted  only  at  specific  times  and  the  model  physics  control  the 
simulation  between  those  times.  Thus,  the  model  physics  can  overcome  biases  in  the  data  that 
are  not  consistent  with  the  model  physics.  This  can  be  an  asset  or  a  liability  depending  on  the 
quality  of  the  model  physics  and  the  data. 

In  the  Newtonian  relaxation  approach,  additional  terms  are  added  to  the  model  equations 
to  assimilate  the  data.  The  model  simulation  is  then  continually  nudged  towards  the  data.  This 
approach  minimizes  the  shock  to  the  model  that  is  often  generated  by  the  reanalysis  procedure  in 
the  discontinuous  data  assimilation  approach,  and  also  permits  data  from  any  time  period  to  be 
easily  used  in  the  data  assimilation  process.  However,  there  are  problems  with  the  Newtonian 
relaxation  approach.  The  presence  of  an  additional  non-physical  term  in  the  model  equations 
means  that  the  result  of  the  simulation  does  not  satisfy  the  true  model  equations  unless  the  data 
assimilation  term  is  zero  (i.e.  the  model  value  and  the  observed  value  of  a  variable  are  identical 
at  a  point).  This  would  not  be  a  significant  issue  if  all  the  data  had  no  errors  and  were 
representative  of  the  model  grid  scale.  Of  course,  that  is  not  the  case.  Therefore,  continuous 
data  assimilation  can  introduce  biases  into  the  model  simulations.  An  example  is  the  previously 
cited  problems  with  rawinsonde  moisture  data.  Continuous  assimilation  of  rawinsonde  moisture 
data  can  introduce  a  dry  bias  into  a  simulation.  Some  modelers  have  addressed  this  problem  by 
using  a  lower  nudging  coefficient  for  the  assimilation  of  rawinsonde  moisture  data  (e.g.  Stauffer 
et  al.,  1991).  Of  course,  this  also  reduces  the  impact  of  the  moisture  data  on  the  entire  simulation. 

Thus,  it  is  not  intuitively  obvious  which  method  would  yield  the  best  climatological 
statistics.  It  is  even  possible  that  one  method  may  produce  better  results  for  some  variables 
while  another  method  will  yield  better  results  for  other  variables.  Experiments  with  both  of  these 
techniques  were  executed  as  part  of  the  current  project.  These  results  provided  guidance  for  the 
selection  pf  the  method  to.be  used  for  the  simulations  over  Korea  and  the  Middle  East. 

The  simulation  strategy  which  was  used  in  the  proof-of-concept  experiments  is 
schematically  illustrated  in  Attachment  2,  This  strategy  represents  one  of  the  options  that  could 
of  been  employed  in  the  current  project.  In  this  approach,  a  coarse  mesh  simulation  (40  km  grid) 
was  initialized  at  0000  UTC  on  the  first  day  of  each  month  for  each  region.  This  simulation  is 
initialized  from  an  optimal  interpolation  (01)  analysis  which  uses  a  first  guess  field  from  a  large 
scale  gridded  analysis  and  all  available  rawinsonde  and  surface  data.  The  large-scale  data  used 
in  the  proof  of  concept  was  the  NCEP  GOI  data.  The  GDAS  Renalysis  data  is  being  used  in  the 
current  project.  The  coarse  mesh  model  was  then  executed  for  24  hours  using  the  large  scale 
interpolated  gridded  analysis  data  (NCEP  GOI  data)  at  12-hour  intervals  to  specify  the  lateral 
boundary  conditions.  The.current  project  is  using  GDAS  Renalysis  data. every  6  hours.  The  data 
from  the  24th  hour  of  this  simulation  was  then  used  as  the  first  guess  for  the  initialization 
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analysis  of  the  next  simulation.  Once  again,  rawinsonde  and  surface  data  were  used  to  update 
the  first  guess.  In  this  scheme,  the  fine  mesh  (10  km  grid)  simulation  was  initialized  at  0000 
UTC  op  the  first  day  of  each  month  by  inteipolating  the  initialization  dataset  of  the  coarse  mesh 
grid.  The  fine  mesh  model  is  then  executed  for  12  hours.  The  lateral  boundary  condition  data 
was  taken  from  the  hourly  output  of  the  40  km  simulation.  The  data  Ifom  the  24th  hour  of  the  1 0 
km  simulation  was  used  to  serve  as  the  first  guess  for  the  initialization  of  the  subsequent  fine 
mesh  simulation.  .All  available  rawinsonde  and  surface  data  are  used  to  update  the  first  guess. 
All  of  the  fields  describing  the  state  of  the  earth's  surface  (soil  moisture,  snow  depth  etc.)  were 
carried  )hrough  the  reinitiahzation  process.  Thus,  the  subsequent  simulation  begins  with  the 
values  existing  at  the  end  of  the  previous  simulation.  If  data  such  as  observed  snow  cover  were 
available,  then  it  probably  would  have  been  beneficial  to  assimilate  it  into  the  simulations. 
Unfortunately  these  types  of  data  are  not  typically  available  in  many  parts  of  the  world. 


The  data  assimilation  scheme  ultimately  selected  for  the  current  project  was  Incremental 
Update  Analysis  (lAU).  MASS  was  used  to  assess  the  lAU  scheme.  The  implementation  of 
lAU  in  mass  required;  (1)  the  modification  of  the  data  preprocessor  so  that  it  can  write  analysis 
increments  rather  than  actual  grid  point  values  when  preparing  for  an  lAU;  and  (2)  modifying  the 
model  so. that  it  can  read  the  LAU  file  written  by  the  preprocessor,  transform  it  to  a  time  tendency 
and  add  this  tendency  to  the  model  variables  at  each  grid  point  during  the  lAU  period.  A  31-day 
series  of  MASS  simulations  extending  fi-om  1  January  to  31  January  1990  was  selected  to 
compare  the  performance  of  the  lAU,  Newtonian  relaxation  and  intermittent  data  assimilation 
(optimal  interpolation)  schemes.  This  time  period  was  chosen  because  MESO  had  already  been 
performed  intermittent  assimilation  and  Newtonian  relaxation  experiments  for  this  period.  A 
schematic  of  the  simulation  strategy  used  in  this  lAU  experiment  is  shown  in  Attachment  3. 
After  the  first  12  hour  simulation  (simulation  A),  each  subsequent  simulation  was  run  for  18 
hours.  The  first  6  hours  of  each  18  hour  simulation  overlapped  with  the. previous  simulation  so 
that  an  lAU  could  be  performed  during  this  period.  The  6-hour  interval  was  considered  to  be  the 
assimilation  period  which  produced  an  initialization  dataset  valid  at  the  end  of  the  interval.  The 
last  12  hours  of  each  18  hour  simulation  (e.g.  simulation  B)  was  evaluated  as  a  forecast.  In  this 
approach,  each  day  consisted  of  two  12-hour  forecasts.  Each  of  these  forecasts  is  initialized  with 
an  lAU  data  assimilation  procedure  with  the  exception  of  the  first  simulation  on  the  first  day 
(simulation  A)  which  is  initialized  with  .an  optimum  interpolation  analysis  of  rawinsonde  data 
employing  an  NCEP  GOI  as  a  first  guess.  The  net  result  was  a  significant  reduction  of  the  biases 
in  the  ipoisture  related  variables  relative  to  the  other  data  assimilation  schemes  while  retaining 
the  observational  control  of  model  drift. 

Attachment  4  shows  all  the  data  assimilation  schemes  considered  for  this  project.  It  was 
determined  that  it  is  best  to  use  a  12  hour  data  assimilation  of  all  available  surface  and 
rawinsonde  data  to  control  model  drift.  To  reduce  the  influences  of  the  biases  in  the  observed 
data  I  Ay  was  implemented  as  the  data  assimilation  scheme.  lAU  is  employed  during  the  first 
three  hours  after  the  assimilation  of  observed  data.  First  lAU  to  calculates  model  and  observed 
data  tendencies.  The  tendencies  are  then  used  to  gently  push  the  simulation  to  the  observed 
solution  without  forcing  the  observed  biases  into  the  model.  The  result  has  been  a  significant 


improvement  in  the  moisture-related  variables  relative  to  the  other  data  assimilation  schemes 
during  the  first  few  hours  of  the  simulations. 

Another  significant  factor  is  the  microphysics  employed.  Our  evaluation  showed  that  to 
optimize  results,  we  needed  to  use  a  complex  prognostic  microphysical  package  which 
incorporates  phase  changes  of  the  moisture  in  the  atmospheric.  Using  simple  traditional 
diagnostic  microphysical  routines  significantly  degraded  the  moisture  related  climate  statistics. 
Unfortunately,  the  more  optimal  climate  run  configurations  are  also  more  computationally 
intensive.  It  was  .jointly  decided  that  the  better  results  justified  the  increase  in  run  time. 

The  final  model  configuration  selected  was: 

(a)  Data  assimilation:  Assimilate  all  available  rawinsonde  and  surface  data  using  lAU 
every  12  hours  using  a  3  hr  assimilation  period 

(b)  Model  physics  for  the  coarse  mesh  (40  km):  The  Blackadar  planetary  boundary  layer 
scheme,  a  modified  version  the  Kou-  Anthes  cumulus  scheme  and  diagnostic  moisture  physics 
were  all  chosen  for  computational  efficiency. 

(c)  Model  physics  for  the  fine  mesh  (10  km):  The  TKE  planetary  boundary  layer  scheme; 
the  Fritsch-Chappell  cumulus  scheme  and  prognostic  moisture  mixed  phase  explicit  physics  with 
prediction  variables  for  cloud  ice,  cloud  water,  rain  water  and  snow  layer  were  all  chosen  to 
optimize  ^he  results  at  the  higher  resolution. 

It  would  be  optimal  to  use  the  TKE  boundary  scheme  for  both  runs.  However,  TKE  takes 
longer  although  it  is  somewhat  better  for  producing  accurate  low-level  fields.  The  choice  of 
Blackadar  is  intended  to  speed  up  the  production  process.  The  choice  of  Blackadar  has  less  of  an 
impact  on. the  40  km  results.than  on  the  10  km  results.  Fritsch-Chappell  handles  the  small-scale 
convective  precipitation  better  than  the  Kuo-Anthes  at  the  10  KM  resolution.  The  explicit 
physics  scheme  is  more  computationally  expensive  but  it  does  produce  superior  results, 
especially  at  the  higher  resolution  of  10  Km. 

An  example  of  the  geographical  domains  for  the  coarse  (40  km)  and  fine  (10  km)  grids 
for  each  region  is  presented  in  Attachments  5  and  6. 
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2.5  CLIMOD  Simulation  Execution  Procedures 

The  general  procedure  for  executing  a  CLIMOD  simulation  begins  with  installing  the 
required  .files  and  directory  structure  and  paths.  This  must  be  done  simultaneously  for  . each 
platform  that  is  being  used.  Currently  five  systems  are  in  use  so  simulations  are  done  for  five 
separate  years.  The  simulations  are  initiated  and  controlled  by  a  main  script,  once  the  setup  is 
completed.  The  only  limitation  to  the  number  of  simulation  days  that  can  be  executed  at  one 
time  is  the  available  disk  space  to  store  the  output  data.  Once  a  simulation  day  is  completed  the 
person  monitoring  the  runs  is  notified 

2-6  Quality  Control  Procedures  During  CLIMOD  Execution 

Extensive  quality  control  procedures  have  been  developed  using  an  HTML-based  system 
on  an  in-house  web  site  (Attachments  7  and  8).  This  allows  continuous  monitoring  of  the 
individual  days  both  across  grid  and  at  point  locations.  Once  a  simulation  day  is  completed  the 
person  monitoring  the.  runs  is  notified  so  they  can  perform  the  quality  control  checks  on  the 
output.  There  are  scripts  that  automate  checking  for  proper  file  size,  missing  output,  and  proper 
format.  The  next  step  is  to  look  at  the  plotted  gridded  data  for  meteorological  soundness. 
Finally,  selected  individual  station  point  locations  of  simulated  observations  are  compared  with 
actual  pbaerved  observations  for  the  same  time  and  point.  If  the  output  data  is  not 
meteorologically  sound  investigation  then  takes  place  to  determine  the  cause.  This  occurs  quite 
often  when  CLIMOD  is  first  applied  to  a  new  location,  as  certain  parameter  values  such  as  initial 
soil  moisture  content  may  have  to  be  modified.  The  modifications  continue  until  the  model  runs 
consistently,  which  has  taken  about  a  week  for  the  two  theaters.  Occasionally  the  run  stops 
because  of  such  things  as  bad  atmospheric  data  or  network  problems.  A  checklist  for 
determining  the  reason  the  run  stopped  and  corrective  action  procedures  have  been  established 
for  each  problem.  A  semi-automated  restarting  procedure  has  been  established  and  is 
continuously  being  improved. 
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3.0-  Computational  Issues 


Several  computational  issues  had  to  be  faced  concerning  both  the  platform  and  the 
operating  system  required  to  meet  a  climate  produaion  capability. 

3.1  Computational  Platform 

Extensive  research  of  the  latest  workstation  technology  was  conducted  to  define  the 
optimum  workstation  and  operating  system  currently  available  within  the  cost  constraints  of  the 
project  budget.  Based  on  cost  considerations,  the  decision  was  to  use  DEC  Alpha  "clones".  The 
specific  systems  used  were  three  500  MHz  Alpha,  and  two  600  MHz  Alpha  processor  systems. 
As  is  typical  when  using  new  systems,  several  hardware  and  software  issues  related  to  system 
compatibility  had  to  be  resolved  in  order  to  use  the  clones  effectively,  but  they  were  resolved. 

3.2.  Operating  System. 

In  order  to  maximize  the  productivity  of  the  resources  available  for  the  project,  MESO 
obtained  and  installed  the  Redhat  Linux  operating  system.  Redhat  Linux  is  considerably  lower 
in  cost 'for  the  same  expected  performance  level  than  a  system  utilizing  operating  systems  from 
"brand"  name  vendors  such  as  Digital  UNIX.  A  Microway  Linux  Fortran  compiler  for 
compiling  MASS  was  also  purchased  and  tested.  This  compiler  is  considerably  less  expensive 
than  the  Digital  Fortran  compiler.  Several  problems  had  to  be  overcome  associated  with  the 
Linux  operating  system.  Even  though  MASS  was  successfully  compiled  and  executed  using  the 
Microway  compiler  under  Linux,  the  performance  of  MASS  never  came  close  to  the 
performance  it  can  achieve  when  compiled  using  the  DEC  Compiler.  MESO  worked  with 
Microwpy  in  an  attempt  to  find  the  solutions  to  optimize  the  performance  of  MASS,  (under  the 
Microway  compiler)  with  limited  success.  The  final  best  solution  at  this  time  was  found  to  be 
compiling  under  the  Digital  operating  system  using  a  DEC  Fortran  compiler  and  then  porting  the 
compiled  MASS  program  to  the  Linux  systems. 


4.0  -Methods  to  Composite  Climate  Data 


It  was  necessary  to  develop  a  post-simulation  method  for  blending  observations  and 
model  generated  data.  The  objective  is  to  develop  a  method  that  will  yield  a  single  seamless  and 
consistent  database  in  such  a  way  that  the  simulated  statistics  surrounding  sites  where  observed 
data  is  available  closely  reflect  the  observed  data.  At  the  same  time,  the  statistics  in  regions 
without  observed  data  are  specified  from  the  model-simulated  data.  This  may  appear  to  be 
somew|iat  redundant  since  the  intention  is  to  assimilate  all  or  most  of  the  available.observed  Hata 
into  the  model  simulations.  Thus,  the  model-generated  database  will  already  be  a  type  of  blend 
betwee^  the  observational  data  and  the  model  physics.  However,  the  model  is  capable  of 
simulating  only  a  portion  of  the  variance  of  each  parameter  because  it  cannot  simulate 
atmospheric  features  that  are  below  the  resolvable  scale  of  the  model  grid.  Of  course,  using  a 
grid  mesh  finer  than  10  km  could  eliminate  this  problem.  However,  that  is  not  computationally 
feasible.  Sub-grid  scale  parameterization  schemes  attempt  to  account  for  some  of  these  features 
but  they  do  not  model  all  of  the  variability  found  at  small  scales.  This  is  particularly  true  for 
variable^  such  as  precipitation,  which  exhibit  a  great  amount  of  small-scale  variability. 
Importantly,  many  extreme  events  (e.g.  isolated  heavy  rains  associated  with  thunderstorms)  are 
associated  with  small-scale  variability.  The  puipose  of  the  post-simulation  blending  technique 
is  to  create  a  single  database  that  will  contain  most  or  all  the  attributes  of  the  observational 
statistics  where  they  are  available  and  blend  smoothly  to  a  model-generated  database  in  other 
areas. 


The  formulation  of  the  blending  method  is  based  on  the  principle  that  the  statistical 
database  should  be  heavily  influenced  by  the  observational  data  near  the  actual  observing  sites 
and  blend  smoothly  into  the  model  statistics  away  from  observational  data  sites.  The  final 
compos^  climatology  should  also  contain  the  full  distribution  of  data,  not  just  means  and 
variances,  so  that  a  large  suite  of  statistical  parameters  can  be  obtained. 

The  blending  procedure  uses  a  weighted  average  of  model  and  observational  data.  The 
weights  depend  upon  the  spatial  statistics  of  each  variable.  In  some  regions  and  in  some 
directions  the  extent  of  the  region  of  high  spatial  correlation  will  be  larger.  In  these  cases  the 
impact  of  the  observational  data  will  likewise  reach  further.  The  spatial  correlation  depends  on 
a  variety  of  factors  including  the  terrain,  the  atmospheric  dynamics,  and  the  time  of  year.  An 
estimate  of  the  spatial  correlation  is  obtained  by  inter-comparing  the  large  grid  simulations,  the 
nested  grid  simulations  and  observational  data  where  available. 

4.1  Results  of  Method  Investigation 

Results  from  Saint  Louis  University's  (SLU)  evaluation  determined  there  are  four  distinct 
method^  to  blend  or  composite  the  observed  data.  All  are  based  fundamentally  on  some  form  of 
distant  weighting  from  the  observed  site.  The  difference  in  the  methods  is  in  the  way  that  the 
autocorrelation  function  relating  the  distance  and  the  weight  is  determined.  The  four  methods 
are; 
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(a)  Imposed  Autocorrelation.  The  distance  autocorrelation  is  calculated 
based  on  derived  relationships  and  reasonable  meteorological  assumptions  associated  with  such 
factors  as  terrain,  slope  and  proximity  to  water. 

(b)  Observation  Derived  Autocorrelation:  The  autocorrelation  function  is 
calculated  based  upon  the  relationship  found  between  observed  sites.  The  major  disadvantage  is 
the  lack  of  a  sufficient  density  of  sites. 

(c)  Model  Derived  Autocoirelation:  The  autocorrelation  function  is 
calculated  based  upon  the  relationship  found  between  model  grid  points.  The  advantage  is  there 
is  adequate  grid  coverage.  The  disadvantage  is  the  method  will  increase  computational 
requirements. 


(d)  A  Composite  Autocorrelation:  This  would  use  a  blend  of  the  three 
methodjs.  It  is  likely  that  this  would  produce  the  best  method,  but  will  require  more  research 
beyond  the  scope  of  this  project. 


The  decision  was  made  to  proceed  with  the  imposed  autocorrelation  method  for  several 
reasons.  The  primary  reason  it  is  that  any  reasonable  autocorrelation  method  is  likely  to  produce 
acceptable  results.  Secondly,  the  method  is  less  computationally  intensive  than  the  other 
method?.  Finally,  with  a  "better"  method,  the  expected  improvement  would  be  slight  and  not 
worth  the  computational  increases.  However,  this  is  an  area  where  further  improvements  could 
probably  ^e  made  to  the  method  with  additional  research. 

Analysis  of  model  results  indicates  that  the  model  statistics  differ  by  only  a  small  amount 
so  that  fhe  blending  should  only  make  small  adjustments.  However,  quantifying  the  differences 
in  blending  algorithms  is  an  area  where  further  improvements  could  probably  be  made  to  the 
method  with  additional  research.  As  part  of  the  current  jjroject  tasks,  the  estimated  improvement 
in  the  climatological  statistics  that  can  be  obtained  through  the  use  of  the  blending  approach  will 
be  rigorously  tested. 

4.2  Description  of  the  Method. 

The  first  step  in  formulating  the  blending  algorithm  is  to  determine  the  weights  for  the 
model  ^nd  blended  observational  data.  The  blended  value  is  related  to  the  model  and 
observations  by: 


(4.1) 

Where  ZJ  represents  the  variable  be  blended  at  the  grid  point  i  (in  this  case  temperature,  but  any 
variable  could  be  substituted),  w  is  the  weight  given  to  the  model  and  observations  at  each 
respective  gridpoint,  and  the  superscripts  B,M,  and  Obs  denote  the  blended,  model  and 
observational  data  respectively.  The  subscript  i  refers  to  the  grid  point  and  the  subscript  j  refers 


21 


to  the  station  location.  The  weights  for  the  observational  data  are  obtained  by  solving  the 
standard  system  of  equations  associated  with  particular  algorithm. 

The  version  considered  initially  was  one  where  the  weights  were  determined  by  the 
variation  in  the  terrain  from  the  observing  site  to  the  center  of  the  grid  point  under  consideration. 
In  particular,  the  weights  are  proportional  to: 

1 


(4.2) 


Where  D  js  the  distance  from  the  observation  point  defined  as: 


D  =  [(n,~n(l)fdl 


(4,3) 


nij  is  the  normal  vector  to  the  surface  at  the  observation  site  and  n(])  is  the  normal  to  the  surface 
along  the  great  circle  path  from  the  observing  site  to  the  model  grid  point. 


The  weight  for  the  model  value  (at  the  grid  4)oint)  is  determined  by  the  condition: 


W  ■  +  ^  w..  =  1 

mi  /  j  ij  ^ 


(4.4) 


There  are  several  degrees  of  freedom  in  this  procedure,  which  include  the  overall 
proportionality  constant  for  Wy  and  the  resolution  at  which  the  normal  vector  integration  is 
determin^dL  JBoth  of  these  effects  act  to  “control”  the  reach  of  the  observational,  data. 

The  second  method  used  to  blend  the  observations  was  a  hybrid  of  the  above  scheme.  In 
this  scheme  the  weights  Wy  are  modified  by  the  weights  used  in  the  MASS  model  initialization 
scheme  to  include  observations.  This  optimal  interpolation  scheme  determines  weights  from  the 
following  system  of  equations: 


(4.5) 


Where  P(ry)  is  the  assumed  correlation  between  site  i  and  site  j  and  P(roj)  is  the 
correlation  between  site  j  and  the  model  grid  point.  The  correlation  is  distance  based  and  is 
given  by: 


P{r)  =  cr[(l  +  s)e 


(4.6) 

Where  r  is  horizontal  distance  between  the  two  observations  to  be  correlated,  R  is  the  radius  of 
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influence  and  o,  s,  and  A.  are  adjustable  parameters.  More  details  are  given  in  the  MASS  model 
reference  guide  (1994).  These  weights  are  then  multiplied  to  the  terrain-based  weights  to  obtain 
the  final  blending  weights. 

Both  of  these  blending  algorithms  were  applied  to  one  month  of  January  simulations  for 
upstate  New  York.  The  differences  in  the  resulting  statistics  were  quite  small  and  the  differences 
between  the  blended  and  original  model  statistics  were  also  quite  small.  Most  of  the  changes 
appeared  to  be  bias  corrections  in  the  neighborhood  of  observational  data. 

The  above  blending  algorithms  do  not  strictly  follow  the  Optimal  Interpolation  (01) 
procedqres.  The  difference  between  the  two  algorithms  is  in  how  the  model  data  is  incoiq^orated 
into  the  blending.  Above  the  model  was  treated  as  just  another  observation,  but  in  01  the  model 
is  treated  as  a  first  guest  (or  background  estimate).  It  is  possible  that  an  01  scheme  will  produce 
better  results.  The  described  algorithms  are  being  adapted  to  follow  an  01  treatment  of  the 
problerp  for  comparison  with  results  obtained  from  the  described  method. 
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5.0-  Determining  Simulated  Climate  Data  Quality 

The  goal  of  this  task  is  to  develop  an  automated  procedure  to  objectively  measure  the 
quality  of  model-simulated  (CLIMOD)  climatologies  in  both  data  rich  and  data  sparse  regions 
and  thus  provide  an  indication  of  the  reliability  of  the  CLIMOD  statistics.  The  CLIMOD 
technique  has  the  capability  to  generate  climate  statistics  anywhere  in  the  world.  However,  the 
statistics  will  not  be  of  equal  quality  due  to  many  factors  including  the  amount  and  t3q)e  of  input 
data  av^able  for  the  numerical  simulations ,  and  the.  particular  model  physics  and  resolution 
(which  may  vary  due  to  availability  of  computational  resources)  used  to  generate  the 
climatological  dataset. 

As  SLU  has  progressed  they  found  that  there  are  two  separate  subdivisions  or  subtasks  to 
this  effort.  The  first  subtask  is  the  assessment  or  validation  of  the  climate  statistics  in  the  vicinity 
of  observations.  At  these  locations  a  direct  measure  of  the  accuracy  of  the  statistics  is  available 
and  hence  the  quality  can  be  determined.  The  second  subtask  is  an  estimate  of  the  confidence  of 
the  clirhate  statistics  at  the  locations  where  observational  data  are  not  available.  Confidence  is 
NOT  a  direct  measure  of  accuracy  but  a  indirect  measure  of  the  expected  quality  of  the  data.  So, 
for  the  very  data  sparse  region  the  confidence  may  be  indicated  as  low  only  because  we  have  no 
objective  way  for  determining  the  quality  of  the  simulated  statistics.  It  is  possible  that  the 
statistics  in  a  low  confidence  location  are  actually  very  good,  there  is  just  know  way  of 
objectively  determining  the  representativeness  of  the  statistics 

5. 1  Assessment  of  the  Method 

The  assessment  procedure  is  based  on  the  procedure  developed  in  the  proof-of-concept 
phase  of  this  research.  The  assessment  must  be  done  against  independent  data  so  the  procedure 
has  two  steps.  The  first  step  is  to  identify  and  withhold  certain  observational  data  from  the 
simulation  procedure  and  the  second  is  the  comparison  of  these  withheld  data  to  the  model 
statistics. 

To  identify  which  stations  are  suitable  for  the  assessment,  SLU  developed  a  clustering 
algorithm  that  groups  stations  of  similar  characteristics.  From  each  of  the  clusters,  one  or  more 
stations  are  chosen  to  be  withheld  based  on  geographical  and  length/quality  of  the  record.  For 
practical  purposes  a  minimum  of  10  and  not  more  than  30  stations  should  be  withheld. 

The  assessment  procedure  established  in  the  proof-of-concept  was  designed  to  establish  the 
accuracy  of  the  model  statistics  at  the  model  resolution  fi'om  station  observations.  The  technique 
identifies  three  factors  which  can  be  used  to  define  the  error  in  the  climate  statistics;  (1)  the  bias, 
(2)  the  difference  in  standard  deviations  between  model  and  observations,  and  (3)  the  difference 
in  shapes  between  the  model  and  observational  distributions. 

These  tte  can  be  combined  into  a  single  quantity  that  provides  an  easy  assessment  of  the 
climate  statistics.  While  this  assessment  procedure  has  been  carried  out  for  the  upstate  New 
York  study,  the  Korea  or  Mideast  domains  have  yet  to  be  done.  Examples  of  these  quantities  are 
found  in  the  technical  reports  for  the  upstate  New  York  January  and  July  simulations. 
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5.2  Confidence  Estimation  Method 


The  estimation  of  the  confidence/quality  of  the  model-generated  statistics  in  data  void 
regions  requires  alternate  approaches  based  on  the  dynamics  of  the  region.  In  these  regions,  use 
of  satellite  data,  error  models,  analogs,  and  interpolation  of  known  errors  are  being  examined. 
All  of  these  approaches  have  weaknesses.  Satellite  information  is  often  large  scale  so  much  of 
the  mesoscale  variability  is  lost.  Consequently,  satellite  data  can  only  be  used  to  provide 
consistency  checks.  Error  models  rely  on  statistically  and  empirically  determined  parameters 
obtained  from  the  numerical  model.  Thus,  error  models  provide  no  external  validation  of  data. 
Analogs  assume  that  accuracy  in  one  region  translates  into  accuracy  in  another  region  with 
similar  characteristics.  Interpolation  of  known  errors  (from  the  observational  sites)  assumes  that 
the  error  structure  function  is  relatively  large  scale.  The  upstate  New  York  results  suggest  the 
error  structure  is  mainly  large  scale  but  further  investigations  are  planned  for  the  Korean 
simulations. 

Determining  the  confidence  of  the  climate  statistics  in  data  sparse  regions  may  include 
several  of  these  procedures  and  possibly  others.  Consequently,  the  final  confidence  product 
generated  should  be  flexible  enough  to  include  a  variety  of  procedures. 

6.0  -Analysis  of  CLIMQD  Quipui 

The  global  grid  point  analysis  data  from  the  NCAR/NCEP  Global  Data  Assimilation 
System  (GD AS)  Reanalysis  project  and  all  available  rawinsonde  and  surface  data  were  used  as 
prescribed  in  the  procedures  outlined  in  Section  2.1  to  create  a  3200  km  by  3200  km  40  km 
resolutipn  and  a  800  km  by  800  km  10  km  resolution  hourly  3-dimensional  dataset  for  the  basic 
model  variables  (pressure,  temperature,  moisture,  wind)  for  alO-year  period  (1973-82)  for  the 
following  2  regions  of  interest;  Korea  and  the  Middle  East.  The  hourly  and  climate  data  is  being 
archived  for  each  region  on  Exabyte  high  density  8  mm  tape.  The  process  of  evaluating  the 
output  is  still  in  progress  at  this  time.  The  results  of  the  climate  statistics  shown  here  should  be 
considered  as  early  results  and  not  conclusive.  The  quality  of  the  10-year  datasets  are  being 
objectively  and  subjectively  evaluated  in  each  region  using  all  of  the  available  observational 
data. 


6.1:  Evaluation  of  Korean  Simulated  Climate  Data 

The  results  for  Korea  are  represented  by  Osan  Air  Base  (RKSO)  located  at  37.08  N  and 
127.03  E  with  an  elevation  of  12  meters. 

a.  Pressurp 

The  proof-of-concept  study  conducted  on  the  eastern  Great  Lakes  region  showed  two 
problems  with  the  modeled  pressures.  One  problem  was  an  approximately  4  mb  bias  in  the 
model  pressures  (both  for  10  km  and  40  km  models).  This  consistent  bias  was  due  to  imbalances 
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introduced  by  the  lateral  boundary  conditions.  Consequently,  a  mass  imbalance  was  established 
with  more  mass  leaving  the  model  domain  than  entered  it.  This  caused  the  modeled  pressure  to 
have  a  low  bias.  The  lateral  boundary  condition  algorithms  were  adjusted  to  ensure  the  lateral 
boundary  conditions  are  in  mass  balance,  which  helped  to  eliminate  this  bias. 

The  second  problem  highlighted  in  the  proof-of-concept  project  was  an  absence  of  the 
diurnal  pressure  cycle.  The  reason  for  this  was  the  fact  that  the  diurnal  solar  thermal  tides  have 
never  bpn  built  into  mesoscale  models.  A  diurnal  pressure  correction  was  added  as  a  post 
processing  step  to  the  pressure  field.  The  results  to  this  point  indicate  that  there  is  considerable 
improvement  in  the  model  pressure.  The  results  in  Figure  1  are  for  Osan  Air  Base.  The  figure 
depicts  that  the  MASS  model  10  Km  modeled  pressure  curve  clearly  captures  the  diurnal 
pressure  cycle. 


—  SLP10KJVI  -^SLPOBS  — SLP40KIVI 


SLP  (mb)  RKSO  Jan  73 


Figure  .1.  Hourly  mean  sea  level  pressure  (SLP)  comparison  for  the  observed  (OBS),  10  Km 
and  40  Km  diurnal  pressure  curve  for  the  month  of  Januaiy  1973  at  Osan  Air  Base  (PKSO). 

b.  Temperature. 

The  temperature  comparison  was  consistent  with  the  results  from  the  proof-of-concept 
research  as  it  again  proved  to  be  the  best-modeled  element.  One  difference  noted  was  that  the  10 
Km  results  were  improved  over  the  proof-of-concept  results.  The  reason  for  this  was  a  change  in 
data  assimilation  strategy  that  included  using  rawinsonde  data  for  the  10  Km  as  well  as  40  Km 
domains.  Also,  the  PBL  scheme  was  change  to  TKE  for  the  10  Km  domain.  The  combination  of 
change^  to  the  PBL  scheme  and  data  assimilation  resulted  in  very  representative  temperature 
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TEMP  (F)  RKSO  Jan  73 


—  T  10KM  — T  OBS  — T  40KIVI 


statistics  for  Osan  for  both  the  10  Km  and  40  Km  domains  as  indicated  in  Figure  2. 

Figure  .2.  Hourly  mean  surface  temperature  (T)  comparison  for  the  observed  (OBS),  10  Km 
and  40  Km  diurnal  pressure  curve  for  the  month  of  January  1973. 

c.  Dewpoint 

The  modeled  dewpoints  are  also  very  representative.  At  40  Km  resolution,  a  warm 
dewpoint  temperature  bias  similar  to  that  found  in  the  proof-of-concept  was  observed.  However, 
the  10  Km  dewpoint  results  were  much  improved  as  indicated  in  Figure  .3.  The  improvement 


-<-DP  10KM  DP  OBS  -^DP  40KM 


DEWPOINT  (F)  RKSO  Jan  73 


Figure  3.  Hourly  surface  dewpoint  (DP)  temperature  comparison  for  the  observed  (OBS),  10 
Km  and  40  Km  diurnal  pressure  curve  for  the  month  of  January  1973. 
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in  the  10  Km  results  is  attributed  to  the  change  in  data  assimilation  strategy  and  using  the  TKE 
PBL  scheme. 

d.  Wind  Speed. 

The  analysis  of  the  wind  speed  at  Osan  in  Figure  4  reveals  the  observed  wind  speeds 
were  slightly  greater  during  daylight  hours  than  the  simulated  wind  speeds.  The  10  Km 
simulation  was  somewhat  better  than  the  40  Km  domain  at  capturing  the  daytime  maximum  that 
occurs  in  the  mid  afternoon.  The  nighttime  simulated  wind  speeds  were  higher  than  observed. 
This  is  likely  due  to  the  inability  of  the  model  to  capture  shallow  inversions  at  night  resulting  in 
higher  modeled  wind  speeds. 


-^SP  10KM  ^SP  OBS  -^SP  40KM 


WIMDSPEED  KTS  RKSO  Jan  73 


Figure.  4.  Hourly  mean  surface  wind  speed  (SP  in  knots)  comparison  for  the  observed  (OBS), 
10  Km  and  40  Km  climate  statistics  for  the  month  of  January  1973  at  Osan  AB  (RKSO). 

e.  Precipitation 

The  frequency  of  precipitation  compared  reasonably  well  with  the  observed  at  both  the 
40  Km  and  the  10  Km  resolution  (Figure  5).  The  frequency  of  precipitation  is  defined  here  as 
the  fraction  of  the  days  in  the  month  for  which  rain  occurred  atthe  indicated  hour  of  the  day.  No 
clear  djurnal  pattern  was  observed  in  either  the  observed  or  the  simulated  data,  and  it  was  not 
clear  that  the  10  Km  modeled  data  was  more  representative  than  the  40  Km  resolution  data. 
Howeyer,  it  is  clear  the  40  Km  resolution  simulations  produced  a  smoother  (less  noisy)  curve 


28 


than  either  the  10  Km  or  observed  data.  It  is  likely  that  all  the  precipitation  curves  will  smooth 
as  additional  years  are  analyzed. 


PCP10KM  — PCPOBS  -^PCP40KW1 


PRECIP  RKSO  Jan  73 


Figure  5.  Frequency  of  precipitation  (PCP)  occurrence  comparison  for  the  observed,  10  Km  and 
40  Km  climate  statistics  for  the  month  of  January  1973  Osan  AB  (RKSO). 

f  Cloud  Coverage 


Initial  steps  of  establishing  methods  to  model  clouds,  tops  and  ceilings,  were  undertaken. 
The  first  step  was  evaluating  the  model  derived  total  cloud  coverage  with  the  total  cloud 
coverage  that  is  observed.  The  simulated  total  cloud  coverage  is  an  integrated  value  for  the 
entire  cplumn  assuming  maximum  overlapping  of  levels.  The  cloud  coverage  at  each  level  is 
obtained  by  calculating  the  cloud  fraction  based  upon  a  simple  cloud  fraction-relative  humidity 
relationship.  The  evaluation  of  cloud  coverage  is  particularly  difficult  because  cloud  coverage 
is  not  always  reported  in  the  observation.  However,  a  limited  comparison  of  the  modeled  versus 
observed  cloud  coverage  for  Osan  AB  (Figure  6)  revealed  encouraging  results.  Currently, 
statistics  are  being  generated  for  ceilings  and  cloud  tops  for  both  Korea  and  the  Mideast.  Once 
compiled,  these  results  will  be  compared  with  available  cloud  observations 
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-^CLDIOKM  —CLDOBS  -^CLD40KM 


CLOUD  FRACTION  RKSO  Jan  73 


Figure  6.  Hourly  mean  total  cloud  coverage  (CLD)  comparison  for  the  observed,  10  Km  and  40 
Km  clipiate  statistics  for  the  month  of  January  1973  at  Osan  AB  (RKSO). 


6.2:  Evaluation  of  Mideast  Simulated  Climate  Data 

A  limited  evaluation  of  the  early  results  from  the  Mideast  simulation  revealed  a  significant 
discrepancy  in  the  temperature  and  dewpoint  as  the  Figures  7  and  8  indicate.  The  results  from  5 
January  1973  from  Baghdad  (32.23  N,  43.15  E,  elevation  34  meters)  are  shown  to  demonstrate 
the  problems  and  the  results  from  modifications  to  the  soil  moisture  to  solve  the  problem.  Both 
Figure  7  and  8  are  initial  results  obtained  before  soil  moisture  content  was  adjusted. 

Figure  7  shows  a  strong  overnight  warm  bias  in  the  modeled  temperature.  The  dewpoint  curve  at 
Baghd^  for  5  January  J973  depicted  in  Figure  8  showed  an  even  greater  problem  with  a  model 
warm  bias.  The  possible  causes  for  this  were  examined.  It  was  determined  that  the  most  likely 
case  wa^  an  unrealistically  high  initial  soil  moisture  content.  Since  the  soil  in  that  region  is  very 
dry,  the  model's  initial  soil  moisture  content  was  changed  to  a  near  zero  value.  This  greatly 
improved  the  results.  The  Baghdad  results  depicting  surface  temperatures  and  dewpoints  for  5 
January  1973  using  the  modified  soil  moisture  content  are  shown  in  Figures  9  and  10.  The 
temperature  results  using  the  modified  soil  moisture  content  decreased  the  overnight  simulated 
low  temperature  values.  The  10  Km  temperature  results  were  improved  more  than  the  40  km 
results.  The  modeled  dewpoint  curve  showed  a  remarkable;  improvement  when  the  soil  moisture 
was  adjusted  as  depicted  in  Figure  10. 
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Figure  7.  The  surface  temperature  comparison  for  the  observed  (OBS),  10  Km  and  40  Km 
domains  for  5  January  1973  at  Baghdad,  Iraq  (ORBB)  before  soil  moisture  was  adjusted. 


Figure  8.  Surface  dewpoint  temperature  comparison  for  the  observed  10  Km  and  40  Km  for  5 
January  1 973  at  Baghdad,  Iraq  (ORBB)  before  soil  moisture  was  adjusted. 
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—  TIOKM  -<^TOBS  -*-T40KM 


iCOR  TEMP  (F)  ORBB  5  Jan  73 


Figure  9.  The  surface  temperature  (T)  comparison  for  the  observed  (OBS)  10  Km  and  40  Km  for 
5  January  1973  at  Baghdad,  Iraq  (ORBB)  ^er  soil  moisture  was  adjusted. 


-»-DP  lOKM  -^DP  OBS  -4-DP  40KM 


DEWPOINT  (F)  ORBB  S  Jan  73 


FigurelO.  The  corrected  surface  dewpoint  (DP)  comparison  for  the  observed  (OBS),  10  km  and 
40  km  for  the  5  January  1973  at  Baghdad,  Iraq  (ORBB)  after  soil  moisture  was  adjusted. 
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6.3:  Construction  of  the  Basic  Climate  Data  Sets 

MESO  did  extensive  consultation  with  AFCCC  and  its  customers  to  determine  the  types 
of  basic  variable  (pressure,  temperature,  moisture,  wind  and  precipitation)  statistics  that  are 
desired  in  the  grid  point  climatologies  for  each  of  the  two  geographical  regions.  A  wide  array  of 
products  have  been  produced  from  four  years  of  the  simulated  climate  data  for  both  Korea  and 
the  Mideast  in  order  to  examine  the  quality  of  the  products.  Currenttly,  MESO,  AFCCC  and 
some  of  its  customers,  and  Saint  Louis  University  are  examining  the  data  for  quality. 
Attachment  9  provides  an  example  of  a  January  temperature  climatology  for  the  Korean  region 
produced  from  four  years  of  simulated  data.  A  primary  and  secondary  (i.e.  backup)  archive  of 
10-year  climatological  statistics  data  and  10  years  of  hourly  simulated  data  for  the  coarse  and 
fine  mesh  grids  is  being  constructed  for  each  of  the  two  regions  in  Exabyte  high-density  8  mm 
tape. 
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7.0  -  Summary 


All  indications  are  that  the  CLEMOD  method  will  provide  high  quality  climate  statistics 
for  both  the  Korean  and  Mideast  regions.  However,  the  results  of  the  climate  statistics  shown 
here  should  be  considered  as  encouraging  early  results  and  not  a  totally  conclusive  evaluation. 
The  quality  of  the  10-year  datasets  will  continue  to  be  objectively  and  subjectively  evaluated  in 
each  region  using  all  of  the  available  observational  data. 

The  primary  task  over  the  next  12  months  will  be  to  finish  evaluating  the  quality  of  the 
10  year  datasets.  In  addition,  specific  data  denial  experiments  will  be  performed  to  help  with 
defining  confidence  levels  and  investigation  will  take  place  to  extend  this  method  to  more 
complex  climate  variables. 


ATTACHMENTS 
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Continue  Model  Rmi 


Attachment  2 
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A  schematic  of  the  data  assimilation  strategy  used  for  the  generation  of  climatological  datasets 
during  th^jjroof-of-concept  research. 
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Attachment  3 
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f  Simulation  C 


Attachment  4 
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Attachment  5 
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Attachment  7 
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Attachment  J8 


High  Resolution  Climate  Project:  Korea 
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Attachment  9 


Mean  Surface  Temperature  (F) 
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Source:  Advanced  Climate  Model  (ACMES) 
All  hOUns.  :^,76, 77,79. 


An  example  of  the  simulated  climatological  gridded  data  produced  by  CLIMOD. 
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